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A general n-state directed 'sandpile' model is introduced. The stationary properties of the n-state 
model are derived for n < oo, and analytical arguments based on a central limit theorem show that 
the model belongs to the universality class of the totally asymmetric Oslo model, with a crossover 
to uncorrelated branching process behavior for small system sizes. Hence, the central limit theorem 
allows us to identify the existence of a large universality class of one-dimensional directed sandpile 
models. 
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Sandpile models have attracted much ana- 
lytical attention in recent years, largely due 
to their application to the development of self- 
organized criticality (SOC)[l-4]. Analytical so- 
lutions, scaling arguments and numerics have 
shown that many models share the same criti- 
cal exponents and scaling functions, leading to 
the notion of universality classes for such sys- 
tems, as in equilibrium systems [2, 5, 6]. In the 
following, we map an n-state directed model to 
a random walker problem and, using a central 
limit theorem for dependent random variables 
[7], derive the conditions for scaling and the 
associated critical exponents. The exponents 
are in exact agreement with those derived for 
the Totally Asymmetric Oslo Model (TAOM) 
[8, 9], which is a special case of this general 
model. For small system sizes, we find that the 
model may exhibit different scaling, which cor- 
responds to an uncorrelated branching process, 
with a crossover characterized by an n depen- 
dent crossover length, £„. 

The model considered is an n-state directed 
'sandpile' model. The system exists on a one- 
dimensional lattice with L sites. Each site, i, is 
in one of n states, Zi £ [0, n — 1], which repre- 
sents the number of particles on site i. 

At the beginning of each time step, we add 
one particle to site i = 1. This site may then 
topple a number of times, each toppling redis- 
tributing a particle to the next site, z\ — ► Z\ — 1 
and Z2 — > Z2 + 1. When site i = 2 receives a par- 
ticle it may also undergo topplings, redistribut- 
ing particles to site i — 3, and so on, with par- 
ticles being passed to sites of increasing i. Note 
that when site i = L topples, the redistributed 
particle will leave the system. When all activ- 
ity ceases, a new time step commences. The 
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avalanche size, s, is defined as the total num- 
ber of topplings during a single time step. The 
only restrictions on the toppling rules are: (i) 
A toppling may never cause z% to become neg- 
ative, (ii) If Zi > n — 1, then site i must topple, 
(iii) Each time site i topples, it redistributes 
one particle to site i + 1 only, (iv) Particles are 
conserved in the bulk, only leaving the system 
at the boundary site i = L and entering when 
particles are added to the boundary site i = 1 
at the beginning of each time step, (v) The 
number of topplings a site undergoes is non- 
deterministic for at least one value of z. This 
final restriction disallows deterministic toppling 
rules, which lead to trivial dynamics. 

For the following general discussion, there is 
no need to specify the toppling rules in any 
more detail. Later, when presenting numerical 
results, we will consider a specific implemen- 
tation of a general class of probabilistic top- 
pling rules that satisfies the 'restrictions' (i)-(v) 
above. 

The quantity of interest in a sandpile which 
has reached a non-equilibrium steady state is 
the avalanche-size probability, P(s; L, 1), which 
is the probability of observing an avalanche of 
size s in a system of size L when one particle 
is added to site i = 1. SOC is associated with 
a time-independent avalanche-size probability 
which obeys simple finite-size scaling 

P(s; L, 1) = as~ Tt S (s/bL D ) fors,L»l, (1) 

where a and b are non-universal constants, r 
and D are universal critical exponents, and S? 
is a universal scaling function. The fcth moment 
of the avalanche-size probability is 

oo 

Qg = X> fe P( S ;£,i) 

s=l 

/oo 
ds as k - T ^(s/bL D ) 

= ab 1+k - T G k {L)L D[ - 1+k - T \ (2) 
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where the sum has been approximated by an in- 
tegral and Gk{L) = f™ bL D du u k ~ T &{u). Pro- 
vided that < Gfc(oo) < oo, we have 

Qg = IW= forI»l, (3) 

where -fk = D(l + k — r) is a universal exponent 
and Tfe = ab 1+k ~ T Gk(L) is a non-universal am- 
plitude which is a constant for L 3> 1. Hence, 
the scaling of the moments with system size 
L is a universal feature which is independent 
of particular details of the dynamics if Eq. (1) 
is valid, the approximation to the integral in 
Eq. (2) does not affect the scaling behavior of 

Q^\, and Gk(L) approaches a non-zero con- 
stant for L — > oo. In the following we shall show 
that under a precise set of conditions, Eq. (3) 
will hold with D = 3/2 and r = 4/3. 

Using a simple extension to the Markov ma- 
trix methods used in Ref. [9], it can be shown 
that if the Markov matrix representing the evo- 
lution operator for the model is regular and 
n < oo, then there is a unique stationary state. 
Note that we will only consider toppling rules 
for which the evolution operator is regular and 
so a unique stationary state exists. Following a 
similar calculation as in [9] , we find that in this 
state, the number of particles, Zi, in each site, 
i, is an independent identically distributed ran- 
dom variable with probability p z [10]. Hence, 
we find that the probability of occurrence of a 
configuration {z^ = {zi,z 2 , . . . ,z L } is 

L 

P{zi}=Y[Pz^ ( 4 ) 

1=1 

where the values of p z depend on the details of 
the toppling rules. This is known as a product 
state and has the property that there arc no 
spatial correlations. However, we shall shortly 
argue that there exists a crossover length, £„, 
such that for system sizes L » £„, there are 
temporal correlations which produce non-trivial 
behavior. 
We define 

oo 

Qgn^E^^'™) ( 5 ) 
s=0 

as the fcth moment of the avalanche-size prob- 
ability, P(s; L, to), for a system of size L which 
has received to particles at site i = 1. 

The first moment is easily derived from the 
fact that, in the stationary state, the aver- 
age number of particles which leave the sys- 
tem through the open boundary must equal the 
number of particles added to the system. Each 
of the to particles topples exactly L times, and 

Q ( L = mL, (6) 



implying 71 = D(2 — r) = 1. 

To derive the scaling of higher moments, we 
introduce P(t, s; 1, L, m) as the joint probabil- 
ity that a system of size 1+L which has received 
m particles at site i = 1 undergoes t topplings 
in the first site and s in the remaining L sites. 
We note that since the model is directed and 
the stationary state is a product state, 

P(t, s; 1, L, to) = P(t; 1, m)P(s; L, t). (7) 

At the beginning of each time step we add one 
particle to site i = 1 and it will topple si times 
with probability P(si;l,l). The second site 
therefore receives s\ particles and as a result 
topples S2 times with probability P(s2;l,si). 
The probability of site 2 toppling S2 times, de- 
noted 02 ($2), i s 

00 

M*2)= E^Mi)^!;!,!) (8) 

si=l 

which follows from Eq. (7). If we define 4>i{x) 
as the probability that site i topples x times, 
with <fii(si) = P(si\ 1, 1), then 

00 

&+i(a;) = J2 P ^ i,y)Mv)- (9) 

!/=l 

This describes a discrete random walker on the 
interval [0,oo]. Since activity stops when one 
of the sites topples zero times, there is an ab- 
sorbing boundary at x — 0. The probability of 
hopping from y to x in a single step is given by 
P(x; l,y). If we denote a particular trajectory 
of a walker x(i), i = Q...L, then the corre- 
sponding avalanche size is 

L 

S = ]>>(*) (10) 

i=l 

with x(0) = 1. This corresponds to the area 
under the first L steps of a random walk with 
an absorbing boundary at x = 0. 

Of course, this is a correlated random walker 
and individual steps arc not independent be- 
cause the probability of hopping a certain dis- 
tance varies depending on where the random 
walker is according to P(s; l,m). 

We now define 

00 

s=0 

which is the width of the probability P(s; l,m) 
around the mean value, to. Using a martingale 
theorem [7] , we can show that iff there exists a 
number < M < 00 such that < Q^ m < M 
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FIG. 1: Area under the random walk with an ab- 
sorbing boundary at x = 0. The solid and dashed 
lines are two different trajectories where the former 
has been absorbed at the boundary at the point i* 
on the i axis where x{i*) = 0. The second walker 
survives until L. The area under the curve for the 
first walker is shaded dark gray and the area for the 
second walker is that of the first plus the light gray 
region. 

for all m, then in the limit L — > oo, the distri- 
bution <f>i (x) will converge for large i to that for 
an equivalent independent random walker [10]. 
Hence, we find that all moments scale, 

QgocL^-^/^ii^-t) (12) 

and we can read off the exponents D = 3/2 and 
r = 4/3. 

Hence, we must find the conditions under 

~ (2) 

which Q\ m is non-zero and finite. Note, again, 
that we are assuming the existence of a unique 
stationary state. Consider a site with z par- 
ticles which has received m particles. After 
s topplings have taken place it will contain 
z' = z + m — s particles. Since both z and 
z' must lie between and n — 1, P(s; 1, to) may 
only be non-zero for to — n + 1 < s < m + n — 1. 
Hence Q^ m < i n — 1) 2 > which is finite for 
n < oo. In order to have Q\ m = 0, there must 
be an m for which only an avalanche of size 
s = to is allowed. However, this is only possible 
if the number of topplings a site undergoes on 
receiving a particle is fully deterministic, which 
are trivial dynamics. Hence, if n < oo and the 
toppling rule leads to non-trivial dynamics we 
have < Q ( ^ m <{n- if and Eq. (12) follows. 

The scaling of the TAOM will only be ob- 
served asymptotically for L — > oo. However, 
we hypothesize the existence of an n dependent 
crossover length, £„, such that TAOM scaling 
is observed for L>( n , with different behavior 
for L <C £ n - To see what happens for system 
sizes L <C £m consider a system when a parti- 
cle is added to site i = 1. If the probability, 
p z , that a site is occupied by z particles has 



support for all z <E [0, n — 1], then, for n» 1, 
it is likely that « z « n - 1. If the num- 
ber of particles on a site z is neither close to 

nor n — 1, the propagating avalanches will 
not be sensitive to the medium within which it 
is propagating and the system will be tempo- 
rally uncorrelated. However, as the avalanche 
propagates through the system, fluctuations in 
the number of topplings increase and as each 
subsequent site is less likely to have uncorre- 
lated topplings they will start to feel the fact 
that n is finite. Hence, for small system sizes, 

1 <C L <C £„, the avalanches are uncorrelated 
and will correspond to an uncorrelated branch- 
ing process with exponents D = 2 and r = 3/2 
[11]. For large system sizes, L » £ n , tem- 
poral correlations will emerge and the system 
will crossover to behavior of the TAOM with 
D = 3/2 and r = 4/3. 

We support the above arguments with nu- 
merical data from the following 'typical' real- 
ization: A site i, < Zi < n — 1, which receives 
a particle will topple once with probability 1/2 
and twice with probability 1/4. If Zi = 0, then 
it cannot topple twice and will topple once with 
probability 1/2. If = n — 1, it will topple once 
or twice, each with probability 1/2. It can be 
shown that the support of p z extends over all 
possible states z G [0,n — 1]. We shall com- 
pare numerical results from this model with ex- 
act results from the corresponding uncorrelated 
branching process, which we can calculate ana- 
lytically [11]. 

Figure 2(a) displays measurements of 
the rescaled second moment, Q^\/L 5 / 2 vs. 
L. From the arguments above, we expect 

Ql!i/ l5/2 t0 scale like Ll/2 for 1 « i « ^ 
and approach a constant for L ^> £„, which is 
supported by the numerics. We do not attempt 
a data collapse because, although it is clear 
that the crossover length £„ is a non-decreasing 
function of n, it has a non-universal functional 
dependence. 

We also consider the moment ratios 

_ (s k )(s) k - 2 _ r fc it- 2 

» = (a2> *-i = TfT- ( 13 ) 

For an avalanche-size probability of the form 
Eq. (1), gk will be universal constants, that is, 
they only depend on r and Sf. Figure 2(b) dis- 
plays 53 vs. system size, L. Since Gu{L) is 
only constant for L ^> 1, the measured will 
only converge toward the universal constant for 
large L. These have values 53 = 9/5 for the 
branching process and we measure (73 w 1.29 
for the TAOM. For 1«L< £„, the measured 
values follow the exact result for the uncorre- 
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(a) 




FIG. 2: Numerical results for a typical realization 
with n = 2,4, 8, 64 and exact results for the un- 
corrected branching process (BP). The errors for 
both graphs were calculated using Efron's Jack- 
knife [12] and are smaller than the symbols, (a) 
The rescaled second moment, Q^\/ L 5 ^ 2 , vs. sys- 
tem size, L. For n = 2,4, 8, 64, the rescaled second 
moment approached a constant for t > 1, while 
for the branching process it increases like L 1 ^ 2 . (b) 
The moment ratio, g$, vs. system size, L. For 
n = 2, 4, 8, 64, the moment ratios follow the branch- 
ing process curve for L <C £, n before approaching 
the TAOM value of g 3 w 1.29 for L > £„. 



lated branching process, with a crossover to the 
TAOM curve for L > £„. 



We have shown that a general n-state di- 
rected sandpilc model of self-organized critical- 
ity belongs to the same universality class as the 
totally asymmetric Oslo model, recently solved 
in Ref. [9]. The precise conditions for this uni- 
versality are that the evolution operator is reg- 
ular, the avalanches are non-deterministic, and 
that n is finite. We have argued that there 
is an n dependent crossover length £„, which 
separates uncorrelated branching process expo- 
nents, D = 2 and r = 3/2, for 1 < L < £„ 
from TAOM exponents, D = 3/2 and r = 4/3 
for L £„. This crossover may be considered 
a consequence of temporal correlations emerg- 
ing in the system, which moves it away from 
the uncorrelated branching process, associated 
with mean-field exponents. The conditions for 
a system to be in this universality class have 
been found using a central limit theorem for de- 
pendent random variables [7]. This is the first 
time a technique of applying a central limit the- 
orem to the discrete model has been used to 
explicitly and precisely identify a universality 
class of non-equilibrium self-organized critical 
systems and we expect much new research will 
follow along these lines. 
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